Indicateur de temps

time.ini = Sys.time()

Chargement des libraries et fonctions

source("../../scriptR/functions.R")
## Le chargement a nécessité le package : ggplot2
## Le chargement a nécessité le package : MASS

Enregistrement des fonctions

print(difftime(Sys.time(), time.ini))
## Time difference of 2.274963 secs

Simulations multiples : populations constantes

setwd("multi_new_methis_pop_size_diff/")

mat_cstt_multi = data.frame.stat(seq_ne = c(100,1000,10000), max_gen = 101, max_simu = 100)

mat_cstt_multi = delta.adm.props(mat_cstt_multi)
mat_cstt_multi = delta.adm.gen.0(mat_cstt_multi)
list_mean = data.frame.mean(mat_cstt_multi, max_gen = 101,
                            seq_ne = c(100,1000,10000),
                            col_cstt = c(1:3), col_mean = c(4:68))
mat_cstt_multi_mean = list_mean[[1]]
mat_cstt_multi_var = list_mean[[2]]
remove(list_mean)

mat_cstt_multi_100 = mat_cstt_multi[which(mat_cstt_multi$Generation == 100),]

mat_cstt_multi_gen0 = mat_cstt_multi[which(mat_cstt_multi$Generation == 0),]
write_cstt_plot("delta.mean.adm.props", mat_cstt_multi_mean,
                "pulse_fondation_multi_simu/delta.mean.adm/Ne_cst/delta_mean_ne_cstt_multi", -0.1, 0.01,
                ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
write_cstt_plot("delta.var.adm.props", mat_cstt_multi_mean,
                "pulse_fondation_multi_simu/delta.var.adm/Ne_cst/delta_var_ne_cstt_multi", -0.2, 0.05)

Statistiques à la fondation en fonction de s1.0

####################
write_one_gen(mat_cstt_multi_100, mat_cstt_multi$s1.0[which(mat_cstt_multi$Generation == 0)],
              mat_cstt_multi_100$delta.mean.adm.props,
              mat_cstt_multi_100$Ne, "s1.0", "delta mean adm props",
              t_file = "pulse_fondation_multi_simu/stat_generation_100/ne_cstt_delta_s10",
              t_plot = "Delta mean adm props en fonction de s1.0")
## Saving 7 x 5 in image

####################
write_one_gen(mat_cstt_multi_gen0, mat_cstt_multi_gen0$s1.0,
              mat_cstt_multi_gen0$delta.mean.adm.props,
              mat_cstt_multi_gen0$Ne, "s1.0", "delta mean adm props",
              t_file = "pulse_fondation_multi_simu/stat_generation_0/ne_cstt_delta_s10",
              t_plot = "Delta mean adm props en fonction de s1.0")
## Saving 7 x 5 in image

####################
write_one_gen(mat_cstt_multi_gen0, mat_cstt_multi_gen0$s1.0,
              mat_cstt_multi_gen0$mean.adm.props,
              mat_cstt_multi_gen0$Ne, "s1.0", "mean adm props",
              t_file = "pulse_fondation_multi_simu/stat_generation_0/ne_cstt_mean_adm_s10",
              t_plot = "Mean adm props en fonction de s1.0")
## Saving 7 x 5 in image

####################
write_one_gen(mat_cstt_multi_gen0, mat_cstt_multi_gen0$s1.0,
              mat_cstt_multi_gen0$Fst.s1.adm,
              mat_cstt_multi_gen0$Ne, "s1.0", "Fst.s1.adm",
              t_file = "pulse_fondation_multi_simu/stat_generation_0/Fst_s1_adm_gen0",
              t_plot = "FST s1-adm en fonction de s1.0")
## Saving 7 x 5 in image

####################
write_one_gen(mat_cstt_multi_gen0, mat_cstt_multi_gen0$s1.0,
              mat_cstt_multi_gen0$mean.F.adm,
              mat_cstt_multi_gen0$Ne, "s1.0", "mean.F.adm",
              t_file = "pulse_fondation_multi_simu/stat_generation_0/F_gen0",
              t_plot = "F moyen en fonction de s1.0")
## Saving 7 x 5 in image

####################
write_one_gen(mat_cstt_multi_gen0, mat_cstt_multi_gen0$s1.0,
              mat_cstt_multi_gen0$mean.het.adm,
              mat_cstt_multi_gen0$Ne, "s1.0", "mean.het.adm",
              t_file = "pulse_fondation_multi_simu/stat_generation_0/He_gen0",
              t_plot = "Hétérozygotie moyenne en fonction de s1,0 à la génération 0",
              add_line = TRUE, line_val = c(0.062,0.049),line_txt = c("s1", "s2"),
              line_txt_pos = c(0.0625, 0.0485),
              line_col = c("#DB0073","#6E0B14"))
## Saving 7 x 5 in image

####################
write_ic_var_plot(mat_cstt_multi_mean, mat_cstt_multi_mean$Generation,
                  mat_cstt_multi_mean$delta.mean.adm.props,
                  mat_cstt_multi_mean$Ne, mat_cstt_multi_mean$Ne,
                  "Delta mean adm props en fonction des générations",
                  -0.1, 0, "Ne", "pulse_fondation_multi_simu/delta.mean.adm/Ne_cst/ne_cstt_mean",
                  -0.12, 0, mat_cstt_multi_var$delta.mean.adm.props,
                  "pulse_fondation_multi_simu/delta.mean.adm/Ne_cst/ne_cstt_mean_IC95",
                  ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

####################
write_ic_var_plot(mat_cstt_multi_mean, mat_cstt_multi_mean$Generation,
                  mat_cstt_multi_mean$delta.adm.gen.0,
                  mat_cstt_multi_mean$Ne, mat_cstt_multi_mean$Ne,
                  "Delta adm gen 0 en fonction des générations",
                  -0.12, 0.03, "Ne", "pulse_fondation_multi_simu/delta.adm0.adm/Ne_cst/delta_adm0_ne_cstt_multi",
                  -0.1, 0.03, mat_cstt_multi_var$delta.adm.gen.0,
                  "pulse_fondation_multi_simu/delta.adm0.adm/Ne_cst/delta_adm0_ne_cstt_multi_IC95",
                  ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")))
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

Simulations multiples : populations croissantes

setwd("multi_new_methis_pop_increase_Nu_var/")

seq_nu_multi = sprintf("%.2f", c(0.01, 0.25, 0.49))
seq_ne_ini_multi = c(100)
seq_ne_fin_multi = c(200, 400, 1000, 10000)
seq_combi_multi = as.data.frame(outer(seq_ne_ini_multi, seq_ne_fin_multi, FUN="paste", sep="-"))
seq_combi_multi = as.data.frame.vector(unlist(seq_combi_multi))[[1]]

mat_inc_multi = data.frame.increase.nu(seq_ne_ini = seq_ne_ini_multi,
                                       seq_combi = seq_combi_multi,
                                       seq_nu = seq_nu_multi, max_simu = 100)

mat_inc_multi = delta.adm.props(mat_inc_multi)
mat_inc_multi = delta.adm.gen.0(mat_inc_multi)
seq_ne_ini_multi = c(100)
seq_ne_fin_multi = c(10000)
seq_combi_multi = as.data.frame(outer(seq_ne_ini_multi, seq_ne_fin_multi, FUN="paste", sep="-"))
seq_combi_multi = as.data.frame.vector(unlist(seq_combi_multi))[[1]]

list_mean = data.frame.mean.u(df_stat = mat_inc_multi,
                              max_gen = 101,
                              seq_combi = seq_combi_multi,
                              seq_u = c(0.01, 0.25, 0.49), 
                              col_mean = c(4:68,70:72),
                              col_cstt = c(1:3,69))
mat_inc_multi_mean = list_mean[[1]]
mat_inc_multi_var = list_mean[[2]]
remove(list_mean)

mat_inc_multi_100 = mat_inc_multi[which(mat_inc_multi$Generation == 100),]
write_increase_plot("delta.mean.adm.props", mat_inc_multi_mean, "pulse_fondation_multi_simu/delta.mean.adm/Ne_inc/delta_mean_pop_inc_nu_va_multi", -0.1, 0.01, ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
write_increase_plot("delta.var.adm.props", mat_inc_multi_mean, "pulse_fondation_multi_simu/delta.var.adm/Ne_inc/delta_var_pop_inc_nu_va_multi", -0.2, 0.05)
write_ic_var_plot(mat_inc_multi_mean, mat_inc_multi_mean$Generation,
                  mat_inc_multi_mean$delta.mean.adm.props,
                  mat_inc_multi_mean$U, mat_inc_multi_mean$U,
                  "Delta mean adm props en fonction des générations",
                  -0.1, 0, "U", "pulse_fondation_multi_simu/delta.mean.adm/Ne_inc/ne_inc_mean",
                  -0.12, 0.03, mat_inc_multi_var$delta.mean.adm.props,
                  "pulse_fondation_multi_simu/delta.mean.adm/Ne_inc/ne_inc_mean_IC95",
                  ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

####################

write_ic_var_plot(mat_inc_multi_mean, mat_inc_multi_mean$Generation,
                  mat_inc_multi_mean$delta.adm.gen.0,
                  mat_inc_multi_mean$U, mat_inc_multi_mean$U,
                  "Delta mean adm props en fonction des générations",
                  -0.1, 0, "U", "pulse_fondation_multi_simu/delta.adm0.adm/Ne_inc/delta_adm0_ne_inc_mean",
                  -0.1, 0.03, mat_inc_multi_var$delta.adm.gen.0,
                  "pulse_fondation_multi_simu/delta.adm0.adm/Ne_inc/delta_adm0_ne_inc_mean_IC95",
                  ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")))
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

####################

write_one_gen(mat_inc_multi_100, mat_inc_multi$s1.0[which(mat_inc_multi$Generation == 0)],
              mat_inc_multi_100$delta.mean.adm.props,
              mat_inc_multi_100$U, "s1.0", "delta mean adm props",
              "pulse_fondation_multi_simu/stat_generation_100/ne_inc_delta_s10",
              "delta mean adm props en fonction de s1.0 à la génération 100")
## Saving 7 x 5 in image

Simulations multiples : bottleneck

setwd("multi_new_methis_bottleneck/")

mat_bot_multi = data.frame.bottleneck(lst_ne0 = c(1000), lst_nef = c(1000),
                                      lst_alpha = c(0.1),
                                      lst_nu = c(0.01, 0.25, 0.49),
                                      lst_bottle = c(20, 50, 80), max_simu = 100)

mat_bot_multi = delta.adm.props(mat_bot_multi)
mat_bot_multi = delta.adm.gen.0(mat_bot_multi)
list_mean = data.frame.mean.bottle(df_stat = mat_bot_multi,
                                            max_gen = 101,
                                            seq_combi = c("1000-1000"),
                                            seq_u = c(0.01, 0.25, 0.49),
                                            seq_alpha = c(0.1),
                                            seq_bottle = c(80),
                                            col_mean = c(4:68,74:76),
                                            col_cstt = c(1:3, 69:73))
mat_bot_multi_mean = list_mean[[1]]
mat_bot_multi_var = list_mean[[2]]
remove(list_mean)

mat_bot_multi_100 = mat_bot_multi[which(mat_bot_multi$Generation == 100),]
write_bottle_plot("delta.mean.adm.props", mat_bot_multi_mean,
                  "pulse_fondation_multi_simu/delta.mean.adm/Ne_bot/delta_mean_bottleneck_multi", -0.1, 0.01,
                  ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
write_bottle_plot("delta.var.adm.props", mat_bot_multi_mean,
                  "pulse_fondation_multi_simu/delta.var.adm/Ne_bot/delta_var_bottleneck_multi", -0.2, 0.05)
write_ic_var_plot(mat_bot_multi_mean, mat_bot_multi_mean$Generation,
                  mat_bot_multi_mean$delta.mean.adm.props,
                  mat_bot_multi_mean$U, mat_bot_multi_mean$U,
                  "Delta mean adm props en fonction des générations",
                  -0.1, 0, "U", "pulse_fondation_multi_simu/delta.mean.adm/Ne_bot/ne_bott_mean",
                  -0.12, 0, mat_bot_multi_var$delta.mean.adm.props,
                  "pulse_fondation_multi_simu/delta.mean.adm/Ne_bot/ne_bott_mean_IC95",
                  ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
## Warning: Removed 1 rows containing missing values (geom_point).
## Warning: Removed 1 row(s) containing missing values (geom_path).
## Saving 7 x 5 in image
## Warning: Removed 1 rows containing missing values (geom_point).

## Warning: Removed 1 row(s) containing missing values (geom_path).
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image
## Warning: Removed 1 rows containing missing values (geom_point).

## Warning: Removed 1 row(s) containing missing values (geom_path).

## Warning: Removed 1 rows containing missing values (geom_point).

## Warning: Removed 1 row(s) containing missing values (geom_path).

####################
write_ic_var_plot(mat_bot_multi_mean, mat_bot_multi_mean$Generation,
                  mat_bot_multi_mean$delta.adm.gen.0,
                  mat_bot_multi_mean$U, mat_bot_multi_mean$U,
                  "Delta mean adm props en fonction des générations",
                  -0.1, 0, "U", "pulse_fondation_multi_simu/delta.adm0.adm/Ne_bot/delta_adm0_ne_bott_mean",
                  -0.1, 0.03, mat_bot_multi_var$delta.adm.gen.0,
                  "pulse_fondation_multi_simu/delta.adm0.adm/Ne_bot/delta_adm0_ne_bott_mean_IC95", bot = c(80),
                  ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")))
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

####################
write_one_gen(mat_bot_multi_100, mat_bot_multi$s1.0[which(mat_bot_multi$Generation == 0)],
              mat_bot_multi_100$delta.mean.adm.props,
              mat_bot_multi_100$U, "s1.0", "delta mean adm props",
              "pulse_fondation_multi_simu/stat_generation_100/ne_bot_delta_s10",
              "delta mean adm props en fonction de s1.0 à la génération 100")
## Saving 7 x 5 in image

Multi simulations

print(difftime(Sys.time(), time.ini))
## Time difference of 4.232389 mins

Multi simulations ns modifiés

setwd("ns_modif/new_methis_pop_size_diff_ns_50/")

mat_cstt_ns_50 = data.frame.stat(seq_ne = c(100,1000,10000), max_gen = 101, max_simu = 100)
mat_cstt_ns_50 = delta.adm.props(mat_cstt_ns_50)

setwd("../new_methis_pop_size_diff_ns_75/")

mat_cstt_ns_75 = data.frame.stat(seq_ne = c(100,1000,10000), max_gen = 101, max_simu = 100)
mat_cstt_ns_75 = delta.adm.props(mat_cstt_ns_75)
list_mean = data.frame.mean(mat_cstt_ns_50, max_gen = 101,
                            seq_ne = c(100,1000,10000),
                            col_cstt = c(1:3), col_mean = c(4:67))
mat_cstt_ns_50_multi_mean = list_mean[[1]]
mat_cstt_ns_50_multi_var = list_mean[[2]]

list_mean = data.frame.mean(mat_cstt_ns_75, max_gen = 101,
                            seq_ne = c(100,1000,10000),
                            col_cstt = c(1:3), col_mean = c(4:67))
mat_cstt_ns_75_multi_mean = list_mean[[1]]
mat_cstt_ns_75_multi_var = list_mean[[2]]
remove(list_mean)
p = plot_stat_gen(mat_cstt_multi_mean, mat_cstt_multi_mean$Generation, 
                  mat_cstt_multi_mean$mean.F.adm, mat_cstt_multi_mean$Ne,
                  mat_cstt_multi_mean$Ne,ligne = T,
                  "F en fonction des générations",
                  legd = TRUE)
print(p)

p = p+geom_errorbar(aes(ymin = mean.F.adm - mat_inc_multi_var$mean.F.adm,
                        ymax = mean.F.adm + mat_inc_multi_var$mean.F.adm),
                    width=0.6)
print(p)

p = plot_stat_gen(mat_cstt_ns_50_multi_mean, mat_cstt_ns_50_multi_mean$Generation, 
                  mat_cstt_ns_50_multi_mean$mean.F.adm, mat_cstt_ns_50_multi_mean$Ne,
                  mat_cstt_ns_50_multi_mean$Ne,ligne = T,
                  "F en fonction des générations",
                  legd = TRUE)
p = p + guides(linetype = FALSE)
p = p + labs(color = "Ne")
print(p)

p = p+geom_errorbar(aes(ymin = mean.F.adm - mat_inc_multi_var$mean.F.adm,
                        ymax = mean.F.adm + mat_inc_multi_var$mean.F.adm),
                    width=0.6)
print(p)

p = plot_stat_gen(mat_cstt_ns_75_multi_mean, mat_cstt_ns_75_multi_mean$Generation, 
                  mat_cstt_ns_75_multi_mean$mean.F.adm, mat_cstt_ns_75_multi_mean$Ne,
                  mat_cstt_ns_75_multi_mean$Ne,ligne = T,
                  "F en fonction des générations",
                  legd = TRUE)
p = p + guides(linetype = FALSE)
p = p + labs(color = "Ne")
print(p)

p = p+geom_errorbar(aes(ymin = mean.F.adm - mat_inc_multi_var$mean.F.adm,
                        ymax = mean.F.adm + mat_inc_multi_var$mean.F.adm),
                    width=0.6)
print(p)

Multi simulations ns modifiés

print(difftime(Sys.time(), time.ini))
## Time difference of 4.927537 mins
cor.test(mat_cstt_multi_gen0$mean.het.adm[which(mat_cstt_multi_gen0$Ne == 100)],mat_cstt_multi_gen0$mean.het.adm[which(mat_cstt_multi_gen0$Ne == 1000)], method ="spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  mat_cstt_multi_gen0$mean.het.adm[which(mat_cstt_multi_gen0$Ne == 100)] and mat_cstt_multi_gen0$mean.het.adm[which(mat_cstt_multi_gen0$Ne == 1000)]
## S = 153852, p-value = 0.447
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## 0.07679568
cor.test(mat_cstt_multi_gen0$Fst.s1.adm[which(mat_cstt_multi_gen0$Ne == 100)],mat_cstt_multi_gen0$Fst.s1.adm[which(mat_cstt_multi_gen0$Ne == 1000)], method ="spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  mat_cstt_multi_gen0$Fst.s1.adm[which(mat_cstt_multi_gen0$Ne == 100)] and mat_cstt_multi_gen0$Fst.s1.adm[which(mat_cstt_multi_gen0$Ne == 1000)]
## S = 155810, p-value = 0.5196
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##       rho 
## 0.0650465
cor.test(mat_cstt_multi_gen0$mean.F.adm[which(mat_cstt_multi_gen0$Ne == 100)],mat_cstt_multi_gen0$mean.F.adm[which(mat_cstt_multi_gen0$Ne == 1000)], method ="spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  mat_cstt_multi_gen0$mean.F.adm[which(mat_cstt_multi_gen0$Ne == 100)] and mat_cstt_multi_gen0$mean.F.adm[which(mat_cstt_multi_gen0$Ne == 1000)]
## S = 161242, p-value = 0.7482
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##        rho 
## 0.03245125
cor.test(mat_cstt_multi_gen0$mean.adm.props[which(mat_cstt_multi_gen0$Ne == 100)],mat_cstt_multi_gen0$mean.adm.props[which(mat_cstt_multi_gen0$Ne == 1000)], method ="spearman")
## 
##  Spearman's rank correlation rho
## 
## data:  mat_cstt_multi_gen0$mean.adm.props[which(mat_cstt_multi_gen0$Ne == 100)] and mat_cstt_multi_gen0$mean.adm.props[which(mat_cstt_multi_gen0$Ne == 1000)]
## S = 154658, p-value = 0.4762
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##       rho 
## 0.0719592

Simulations multiples : s1,0 fixé

setwd("s1.0_modif/multi_new_methis_pop_size_diff/s1.0_0.1/")

mat_cstt_multi_s1.0_0.1 = data.frame.stat(seq_ne = c(100,1000,10000), 
                                          max_gen = 101, max_simu = 100)
mat_cstt_multi_s1.0_0.1 = delta.adm.props(mat_cstt_multi_s1.0_0.1)
mat_cstt_multi_s1.0_0.1 = delta.adm.gen.0(mat_cstt_multi_s1.0_0.1)

setwd("../s1.0_0.5/")

mat_cstt_multi_s1.0_0.5 = data.frame.stat(seq_ne = c(100,1000,10000), 
                                           max_gen = 101, max_simu = 100)
mat_cstt_multi_s1.0_0.5 = delta.adm.props(mat_cstt_multi_s1.0_0.5)
mat_cstt_multi_s1.0_0.5 = delta.adm.gen.0(mat_cstt_multi_s1.0_0.5)

setwd("../s1.0_0.75/")

mat_cstt_multi_s1.0_0.75 = data.frame.stat(seq_ne = c(100,1000,10000), 
                                           max_gen = 101, max_simu = 100)
mat_cstt_multi_s1.0_0.75 = delta.adm.props(mat_cstt_multi_s1.0_0.75)
mat_cstt_multi_s1.0_0.75 = delta.adm.gen.0(mat_cstt_multi_s1.0_0.75)

setwd("../s1.0_0.9/")

mat_cstt_multi_s1.0_0.9 = data.frame.stat(seq_ne = c(100,1000,10000), 
                                           max_gen = 101, max_simu = 100)
mat_cstt_multi_s1.0_0.9 = delta.adm.props(mat_cstt_multi_s1.0_0.9)
mat_cstt_multi_s1.0_0.9 = delta.adm.gen.0(mat_cstt_multi_s1.0_0.9)
list_mean = data.frame.mean(mat_cstt_multi_s1.0_0.1, max_gen = 101,
                            seq_ne = c(100,1000,10000),
                            col_cstt = c(1:3), col_mean = c(4:68))
mat_cstt_multi_mean_s1.0_0.1 = list_mean[[1]]
mat_cstt_multi_var_s1.0_0.1 = list_mean[[2]]
remove(list_mean)

list_mean = data.frame.mean(mat_cstt_multi_s1.0_0.5, max_gen = 101,
                            seq_ne = c(100,1000,10000),
                            col_cstt = c(1:3), col_mean = c(4:68))
mat_cstt_multi_mean_s1.0_0.5 = list_mean[[1]]
mat_cstt_multi_var_s1.0_0.5 = list_mean[[2]]
remove(list_mean)

list_mean = data.frame.mean(mat_cstt_multi_s1.0_0.75, max_gen = 101,
                            seq_ne = c(100,1000,10000),
                            col_cstt = c(1:3), col_mean = c(4:68))
mat_cstt_multi_mean_s1.0_0.75 = list_mean[[1]]
mat_cstt_multi_var_s1.0_0.75 = list_mean[[2]]
remove(list_mean)

list_mean = data.frame.mean(mat_cstt_multi_s1.0_0.9, max_gen = 101,
                            seq_ne = c(100,1000,10000),
                            col_cstt = c(1:3), col_mean = c(4:68))
mat_cstt_multi_mean_s1.0_0.9 = list_mean[[1]]
mat_cstt_multi_var_s1.0_0.9 = list_mean[[2]]
remove(list_mean)
lst_tmp = list(mat_cstt_multi_mean_s1.0_0.1,
               mat_cstt_multi_mean_s1.0_0.5,
               mat_cstt_multi_mean_s1.0_0.75,
               mat_cstt_multi_mean_s1.0_0.9)
lst_var = list(mat_cstt_multi_var_s1.0_0.1,
               mat_cstt_multi_var_s1.0_0.5,
               mat_cstt_multi_var_s1.0_0.75,
               mat_cstt_multi_var_s1.0_0.9)
vec_s10 = c("01", "05", "075", "09")

for (i in c(1:4)) {
  write_ic_var_plot(lst_tmp[[i]],lst_tmp[[i]]$Generation,
                    lst_tmp[[i]]$mean.het.adm,
                    lst_tmp[[i]]$Ne,lst_tmp[[i]]$Ne,
                    "He moyen en fonction des générations", 0, 0.07, "Ne",
                    str_c("pulse_fondation_multi_simu/mean.he/Ne_cst/he_mean_ne_cstt_multi_s10_",vec_s10[i]),
                    0, 0.07, lst_var[[i]]$mean.het.adm,
                    str_c("pulse_fondation_multi_simu/mean.he/Ne_cst/he_mean_ne_cstt_multi_s10_",vec_s10[i],"_ic95"),
                    ylab="He moyenne")
  write_ic_var_plot(lst_tmp[[i]],lst_tmp[[i]]$Generation,
                    lst_tmp[[i]]$mean.F.adm,
                    lst_tmp[[i]]$Ne,lst_tmp[[i]]$Ne,
                    "He moyen en fonction des générations", -0.02, 0.005, "Ne",
                    str_c("pulse_fondation_multi_simu/mean.F/Ne_cst/F_ne_cstt_multi_s10_",vec_s10[i]),
                    -0.02, 0.005, lst_var[[i]]$mean.F.adm,
                    str_c("pulse_fondation_multi_simu/mean.F/Ne_cst/F_ne_cstt_multi_s10_",vec_s10[i],"_ic95"),
                    ylab = "F")
  write_ic_var_plot(lst_tmp[[i]],lst_tmp[[i]]$Generation,
                    lst_tmp[[i]]$f3,
                    lst_tmp[[i]]$Ne,lst_tmp[[i]]$Ne,
                    "He moyen en fonction des générations", -0.01, 1, "Ne",
                    str_c("pulse_fondation_multi_simu/f3/Ne_cst/f3_ne_cstt_multi_s10_",vec_s10[i]),
                    -0.01, 1, lst_var[[i]]$f3,
                    str_c("pulse_fondation_multi_simu/f3/Ne_cst/f3_ne_cstt_multi_s10_",vec_s10[i],"_ic95"),
                    ylab="f3")
  write_ic_var_plot(lst_tmp[[i]],lst_tmp[[i]]$Generation,
                    lst_tmp[[i]]$delta.mean.adm.props,
                    lst_tmp[[i]]$Ne,lst_tmp[[i]]$Ne,
                    "Delta mean adm props en fonction des générations", -0.15, 0.1, "Ne",
                    str_c("pulse_fondation_multi_simu/delta.mean.adm/Ne_cst/delta_mean_ne_cstt_multi_s10_",vec_s10[i]),
                    -0.2, 0.1, lst_var[[i]]$delta.mean.adm.props,
                    str_c("pulse_fondation_multi_simu/delta.mean.adm/Ne_cst/delta_mean_ne_cstt_multi_s10_",vec_s10[i],"_ic95"),
                    ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
  write_ic_var_plot(lst_tmp[[i]],lst_tmp[[i]]$Generation,
                    lst_tmp[[i]]$delta.adm.gen.0,
                    lst_tmp[[i]]$Ne,lst_tmp[[i]]$Ne,
                    "Delta mean adm props en fonction des générations",-0.15, 0.1, "Ne",
                    str_c("pulse_fondation_multi_simu/delta.adm0.adm/Ne_cst/delta_adm_ne_cstt_multi_s10_",vec_s10[i]),
                    -0.15, 0.05, lst_var[[i]]$delta.adm.gen.0,
                    str_c("pulse_fondation_multi_simu/delta.adm0.adm/Ne_cst/delta_adm_ne_cstt_multi_s10_",vec_s10[i],"_ic95"),
                    ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")))
}
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Warning: Removed 117 rows containing missing values (geom_point).
## Warning: Removed 117 row(s) containing missing values (geom_path).
## Saving 7 x 5 in image
## Warning: Removed 117 rows containing missing values (geom_point).

## Warning: Removed 117 row(s) containing missing values (geom_path).
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image
## Warning: Removed 117 rows containing missing values (geom_point).

## Warning: Removed 117 row(s) containing missing values (geom_path).

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

remove(lst_tmp)
remove(lst_var)
remove(vec_s10)
setwd("s1.0_modif/multi_new_methis_pop_increase_Nu_var/s1.0_0.1/")
mat_inc_multi_s1.0_0.1 = data.frame.increase.nu(seq_ne_ini = c("100"),
                                                seq_combi = c("100-10000"),
                                                seq_nu = c(0.01, 0.25, 0.49),
                                                max_simu = 100)
mat_inc_multi_s1.0_0.1 = delta.adm.props(mat_inc_multi_s1.0_0.1)
mat_inc_multi_s1.0_0.1 = delta.adm.gen.0(mat_inc_multi_s1.0_0.1)

setwd("../s1.0_0.5/")
mat_inc_multi_s1.0_0.5 = data.frame.increase.nu(seq_ne_ini = c("100"),
                                                seq_combi = c("100-10000"),
                                                seq_nu = c(0.01, 0.25, 0.49),
                                                max_simu = 100)
mat_inc_multi_s1.0_0.5 = delta.adm.props(mat_inc_multi_s1.0_0.5)
mat_inc_multi_s1.0_0.5 = delta.adm.gen.0(mat_inc_multi_s1.0_0.5)

setwd("../s1.0_0.75/")
mat_inc_multi_s1.0_0.75 = data.frame.increase.nu(seq_ne_ini = c("100"),
                                                 seq_combi = c("100-10000"),
                                                 seq_nu = c(0.01, 0.25, 0.49),
                                                 max_simu = 100)
mat_inc_multi_s1.0_0.75 = delta.adm.props(mat_inc_multi_s1.0_0.75)
mat_inc_multi_s1.0_0.75 = delta.adm.gen.0(mat_inc_multi_s1.0_0.75)

setwd("../s1.0_0.9/")
mat_inc_multi_s1.0_0.9 = data.frame.increase.nu(seq_ne_ini = c("100"),
                                                 seq_combi = c("100-10000"),
                                                 seq_nu = c(0.01, 0.25, 0.49),
                                                 max_simu = 100)
mat_inc_multi_s1.0_0.9 = delta.adm.props(mat_inc_multi_s1.0_0.9)
mat_inc_multi_s1.0_0.9 = delta.adm.gen.0(mat_inc_multi_s1.0_0.9)
list_mean = data.frame.mean.u(df_stat = mat_inc_multi_s1.0_0.1, max_gen = 101,
                              seq_combi = c("100-10000"),
                              seq_u = c(0.01, 0.25, 0.49), 
                              col_mean = c(4:68,70:72),
                              col_cstt = c(1:3,69))
mat_inc_multi_mean_s1.0_0.1 = list_mean[[1]]
mat_inc_multi_var_s1.0_0.1 = list_mean[[2]]
remove(list_mean)

list_mean = data.frame.mean.u(df_stat = mat_inc_multi_s1.0_0.5, max_gen = 101,
                              seq_combi = c("100-10000"),
                              seq_u = c(0.01, 0.25, 0.49), 
                              col_mean = c(4:68,70:72),
                              col_cstt = c(1:3,69))
mat_inc_multi_mean_s1.0_0.5 = list_mean[[1]]
mat_inc_multi_var_s1.0_0.5 = list_mean[[2]]
remove(list_mean)

list_mean = data.frame.mean.u(df_stat = mat_inc_multi_s1.0_0.75, max_gen = 101,
                              seq_combi = c("100-10000"),
                              seq_u = c(0.01, 0.25, 0.49), 
                              col_mean = c(4:68,70:72),
                              col_cstt = c(1:3,69))
mat_inc_multi_mean_s1.0_0.75 = list_mean[[1]]
mat_inc_multi_var_s1.0_0.75 = list_mean[[2]]

list_mean = data.frame.mean.u(df_stat = mat_inc_multi_s1.0_0.9, max_gen = 101,
                              seq_combi = c("100-10000"),
                              seq_u = c(0.01, 0.25, 0.49), 
                              col_mean = c(4:68,70:72),
                              col_cstt = c(1:3,69))
mat_inc_multi_mean_s1.0_0.9 = list_mean[[1]]
mat_inc_multi_var_s1.0_0.9 = list_mean[[2]]
remove(list_mean)
lst_tmp = list(mat_inc_multi_mean_s1.0_0.1,
               mat_inc_multi_mean_s1.0_0.5,
               mat_inc_multi_mean_s1.0_0.75,
               mat_inc_multi_mean_s1.0_0.9)
lst_var = list(mat_inc_multi_var_s1.0_0.1,
               mat_inc_multi_var_s1.0_0.5,
               mat_inc_multi_var_s1.0_0.75,
               mat_inc_multi_var_s1.0_0.9)
vec_s10 = c("01", "05", "075","09")

for (i in c(1:4)) {
  write_increase_plot("delta.mean.adm.props", lst_tmp[[i]],
                      str_c("pulse_fondation_multi_simu/delta.mean.adm/Ne_inc/delta_mean_ne_inc_s10_",vec_s10[i]),
                      -0.15, 0.1, name_color = "U",name_line = "Ne",
                      labcolor = "U", labline = "ne0-nef",
                      ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
  # write_increase_plot("delta.adm.gen.0", lst_tmp[[i]],
  #                     str_c("delta_adm_gen0/delta_adm_ne_inc_multi_s10_",vec_s10[i]),
  #                     -0.15, 0.1, name_color = "U",name_line = "Ne",
  #                     labcolor = "U", labline = "ne0-nef",
  #                     ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")))
  write_ic_var_plot(lst_tmp[[i]],lst_tmp[[i]]$Generation,
                    lst_tmp[[i]]$delta.adm.gen.0,
                    lst_tmp[[i]]$U,lst_tmp[[i]]$U,
                    "Delta mean adm props en fonction des générations",-0.15, 0.1, "Ne",
                    str_c("pulse_fondation_multi_simu/delta.adm0.adm/Ne_inc/delta_adm_ne_inc_multi_s10_",vec_s10[i]),
                    -0.15, 0.05, lst_var[[i]]$delta.adm.gen.0,
                    str_c("pulse_fondation_multi_simu/delta.adm0.adm/Ne_inc/delta_adm_ne_inc_multi_s10_",vec_s10[i],"_ic95"),
                    ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")))
}
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

remove(lst_tmp)
remove(vec_s10)
remove(lst_var)
setwd("s1.0_modif/multi_new_methis_bottleneck/s1.0_0.1/")
mat_bot_multi_s1.0_0.1 = data.frame.bottleneck(lst_ne0 = c(1000), lst_nef = c(1000), 
                                               lst_alpha = c(0.1),
                                               lst_nu = c(0.01, 0.25, 0.49),
                                               lst_bottle = c(80), max_simu = 99)
mat_bot_multi_s1.0_0.1 = delta.adm.props(mat_bot_multi_s1.0_0.1)
mat_bot_multi_s1.0_0.1 = delta.adm.gen.0(mat_bot_multi_s1.0_0.1)

setwd("../s1.0_0.5/")
mat_bot_multi_s1.0_0.5 = data.frame.bottleneck(lst_ne0 = c(1000), lst_nef = c(1000),
                                      lst_alpha = c(0.1),
                                      lst_nu = c(0.01, 0.25, 0.49),
                                      lst_bottle = c(80), max_simu = 99)
mat_bot_multi_s1.0_0.5 = delta.adm.props(mat_bot_multi_s1.0_0.5)
mat_bot_multi_s1.0_0.5 = delta.adm.gen.0(mat_bot_multi_s1.0_0.5)

setwd("../s1.0_0.75/")
mat_bot_multi_s1.0_0.75 = data.frame.bottleneck(lst_ne0 = c(1000), lst_nef = c(1000),
                                      lst_alpha = c(0.1),
                                      lst_nu = c(0.01, 0.25, 0.49),
                                      lst_bottle = c(80), max_simu = 99)
mat_bot_multi_s1.0_0.75 = delta.adm.props(mat_bot_multi_s1.0_0.75)
mat_bot_multi_s1.0_0.75 = delta.adm.gen.0(mat_bot_multi_s1.0_0.75)

setwd("../s1.0_0.9/")
mat_bot_multi_s1.0_0.9 = data.frame.bottleneck(lst_ne0 = c(1000), lst_nef = c(1000),
                                      lst_alpha = c(0.1),
                                      lst_nu = c(0.01, 0.25, 0.49),
                                      lst_bottle = c(80), max_simu = 99)
mat_bot_multi_s1.0_0.9 = delta.adm.props(mat_bot_multi_s1.0_0.9)
mat_bot_multi_s1.0_0.9 = delta.adm.gen.0(mat_bot_multi_s1.0_0.9)
list_mean = data.frame.mean.bottle(df_stat = mat_bot_multi_s1.0_0.1,
                                            max_gen = 101,
                                            seq_combi = c("1000-1000"),
                                            seq_u = c(0.01, 0.25, 0.49),
                                            seq_alpha = c(0.1), seq_bottle = c(80),
                                            col_mean = c(4:68,74:76),
                                            col_cstt = c(1:3, 69:73))
mat_bot_multi_mean_s1.0_0.1 = list_mean[[1]]
mat_bot_multi_var_s1.0_0.1 = list_mean[[2]]

list_mean = data.frame.mean.bottle(df_stat = mat_bot_multi_s1.0_0.5,
                                            max_gen = 101,
                                            seq_combi = c("1000-1000"),
                                            seq_u = c(0.01, 0.25, 0.49),
                                            seq_alpha = c(0.1), seq_bottle = c(80),
                                            col_mean = c(4:68,74:76),
                                            col_cstt = c(1:3, 69:73))
mat_bot_multi_mean_s1.0_0.5 = list_mean[[1]]
mat_bot_multi_var_s1.0_0.5 = list_mean[[2]]

list_mean = data.frame.mean.bottle(df_stat = mat_bot_multi_s1.0_0.75,
                                            max_gen = 101,
                                            seq_combi = c("1000-1000"),
                                            seq_u = c(0.01, 0.25, 0.49),
                                            seq_alpha = c(0.1), seq_bottle = c(80),
                                            col_mean = c(4:68,74:76),
                                            col_cstt = c(1:3, 69:73))
mat_bot_multi_mean_s1.0_0.75 = list_mean[[1]]
mat_bot_multi_var_s1.0_0.75 = list_mean[[2]]

list_mean = data.frame.mean.bottle(df_stat = mat_bot_multi_s1.0_0.9,
                                            max_gen = 101,
                                            seq_combi = c("1000-1000"),
                                            seq_u = c(0.01, 0.25, 0.49),
                                            seq_alpha = c(0.1), seq_bottle = c(80),
                                            col_mean = c(4:68,74:76),
                                            col_cstt = c(1:3, 69:73))
mat_bot_multi_mean_s1.0_0.9 = list_mean[[1]]
mat_bot_multi_var_s1.0_0.9 = list_mean[[2]]
remove(list_mean)
lst_tmp = list(mat_bot_multi_mean_s1.0_0.1,
               mat_bot_multi_mean_s1.0_0.5,
               mat_bot_multi_mean_s1.0_0.75,
               mat_bot_multi_mean_s1.0_0.9)
lst_var = list(mat_bot_multi_var_s1.0_0.1,
               mat_bot_multi_var_s1.0_0.5,
               mat_bot_multi_var_s1.0_0.75,
               mat_bot_multi_var_s1.0_0.9)
vec_s10 = c("01", "05", "075","09")

for (i in c(1:4)) {
  write_bottle_plot("delta.mean.adm.props", lst_tmp[[i]],
                    str_c("pulse_fondation_multi_simu/delta.mean.adm/Ne_bot/delta_mean_bottleneck_multi_s10_",vec_s10[i]),
                    -0.1, 0.1,time_bott = c(81), name_color = "U", labcol = "U",
                    ylab=expression(paste(δ[µ["observée"]-µ["théorique"]], sep = "")))
  # write_bottle_plot("delta.adm.gen.0", lst_tmp[[i]],
  #                   str_c("delta_adm_gen0/delta_adm_bott_multi_s10_",vec_s10[i]),
  #                   -0.15, 0.1,time_bott = c(81), name_color = "U", labcol = "U",
  #                   ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")))
  write_ic_var_plot(lst_tmp[[i]],lst_tmp[[i]]$Generation,
                    lst_tmp[[i]]$delta.adm.gen.0,
                    lst_tmp[[i]]$U,lst_tmp[[i]]$U,
                    "Delta mean adm props en fonction des générations",-0.15, 0.1, "Ne",
                    str_c("pulse_fondation_multi_simu/delta.adm0.adm/Ne_bot/delta_adm_ne_bot_multi_s10_",vec_s10[i]),
                    -0.15, 0.05, lst_var[[i]]$delta.adm.gen.0,
                    str_c("pulse_fondation_multi_simu/delta.adm0.adm/Ne_bot/delta_adm_ne_bot_multi_s10_",vec_s10[i],"_ic95"),
                    ylab=expression(paste(δ[µ["g"]-µ[0]], sep = "")),line_bot = TRUE)

}
## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image
## Warning: Removed 35 rows containing missing values (geom_point).
## Warning: Removed 35 row(s) containing missing values (geom_path).

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image
## Warning: Removed 60 rows containing missing values (geom_point).
## Warning: Removed 60 row(s) containing missing values (geom_path).

## Saving 7 x 5 in image
## Scale for 'y' is already present. Adding another scale for 'y', which will
## replace the existing scale.
## Saving 7 x 5 in image

remove(lst_tmp)
remove(vec_s10)
print(difftime(Sys.time(), time.ini))
## Time difference of 9.914185 mins